library(Biobase)
## Loading required package: BiocGenerics
## 
## Attaching package: 'BiocGenerics'
## The following objects are masked from 'package:stats':
## 
##     IQR, mad, sd, var, xtabs
## The following objects are masked from 'package:base':
## 
##     anyDuplicated, aperm, append, as.data.frame, basename, cbind,
##     colnames, dirname, do.call, duplicated, eval, evalq, Filter, Find,
##     get, grep, grepl, intersect, is.unsorted, lapply, Map, mapply,
##     match, mget, order, paste, pmax, pmax.int, pmin, pmin.int,
##     Position, rank, rbind, Reduce, rownames, sapply, setdiff, sort,
##     table, tapply, union, unique, unsplit, which.max, which.min
## Welcome to Bioconductor
## 
##     Vignettes contain introductory material; view with
##     'browseVignettes()'. To cite Bioconductor, see
##     'citation("Biobase")', and for packages 'citation("pkgname")'.
library(survminer)
## Loading required package: ggplot2
## Loading required package: ggpubr
library(survival)
## 
## Attaching package: 'survival'
## The following object is masked from 'package:survminer':
## 
##     myeloma
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✔ tibble  3.2.1     ✔ dplyr   1.1.2
## ✔ tidyr   1.3.0     ✔ stringr 1.5.0
## ✔ readr   2.1.4     ✔ forcats 0.5.1
## ✔ purrr   1.0.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::combine()    masks Biobase::combine(), BiocGenerics::combine()
## ✖ dplyr::filter()     masks stats::filter()
## ✖ dplyr::lag()        masks stats::lag()
## ✖ ggplot2::Position() masks BiocGenerics::Position(), base::Position()
library(DT)
theme_set(theme_gray(base_size = 25))
knitr::opts_chunk$set(dev = "png",
                      dpi = 300,
                      echo = FALSE,
                      fig.height = 6,
                      fig.width = 7, 
                      cache = TRUE)

PATH <- file.path(Sys.getenv("MLAB"), "projects/brcameta/exosome/4t1_brca")
TCGAPATH <- file.path(Sys.getenv("CBM"), "TCGA-GDC/")
DPATH <- file.path(Sys.getenv("CBM"), "otherStudies/RNAseq/2023-03-22-YuhanExosomeBrCa")

do_save <- FALSE

Loading Data

Data Filtering

Censoring: If samples are alive, then they lived at least as long as the day to their last followup. For 5 year threshold, those who died after 5 years, should be right-censored, and labeled to have lived at least as long as ‘days to death’ + 1.

Survfit All subtypes

Adding GSVA data

IR_C_UP

IR_C

IS_C_UP

IS_C

IR_IS_UP

IR_IS

TGFB_C_UP

TGFB_C

Survfit Basal

Adding GSVA data

IR_C_UP

IR_C

IS_C_UP

IS_C

IR_IS_UP

IR_IS

TGFB_C_UP

TGFB_C

Cox Proportional Hazard

TCGA
## Call:
## coxph(formula = Surv(as.numeric(time), vital_status_1) ~ ir_is + 
##     subtype_selected + age_at_index, data = pData(tcga_gsva))
## 
##   n= 1208, number of events= 199 
## 
##                                 coef exp(coef)  se(coef)      z Pr(>|z|)    
## ir_is                      -0.151404  0.859501  0.251669 -0.602  0.54744    
## subtype_selectedBRCA.Basal -0.740883  0.476693  0.283768 -2.611  0.00903 ** 
## subtype_selectedBRCA.Her2  -0.181146  0.834314  0.340328 -0.532  0.59454    
## subtype_selectedBRCA.LumA  -1.035346  0.355104  0.218367 -4.741 2.12e-06 ***
## subtype_selectedBRCA.LumB  -0.610958  0.542830  0.265304 -2.303  0.02129 *  
## age_at_index                0.034946  1.035564  0.005628  6.209 5.33e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                            exp(coef) exp(-coef) lower .95 upper .95
## ir_is                         0.8595     1.1635    0.5248    1.4076
## subtype_selectedBRCA.Basal    0.4767     2.0978    0.2733    0.8313
## subtype_selectedBRCA.Her2     0.8343     1.1986    0.4282    1.6256
## subtype_selectedBRCA.LumA     0.3551     2.8161    0.2315    0.5448
## subtype_selectedBRCA.LumB     0.5428     1.8422    0.3227    0.9130
## age_at_index                  1.0356     0.9657    1.0242    1.0471
## 
## Concordance= 0.697  (se = 0.022 )
## Likelihood ratio test= 67.24  on 6 df,   p=2e-12
## Wald test            = 70.12  on 6 df,   p=4e-13
## Score (logrank) test = 72.23  on 6 df,   p=1e-13
Metabric
## Call:
## coxph(formula = Surv(as.numeric(time), vital_status_1) ~ ir_is + 
##     Pam50_SUBTYPE + AGE_AT_DIAGNOSIS, data = pData(metabric_gsva))
## 
##   n= 1980, number of events= 1143 
## 
##                         coef exp(coef)  se(coef)      z Pr(>|z|)    
## ir_is              -0.140774  0.868686  0.112651 -1.250   0.2114    
## Pam50_SUBTYPEBasal  0.181927  1.199526  0.131293  1.386   0.1659    
## Pam50_SUBTYPEHer2   0.279093  1.321930  0.132840  2.101   0.0356 *  
## Pam50_SUBTYPELumA  -0.216832  0.805065  0.114617 -1.892   0.0585 .  
## Pam50_SUBTYPELumB   0.113128  1.119775  0.117410  0.964   0.3353    
## Pam50_SUBTYPENC     0.106648  1.112543  0.459897  0.232   0.8166    
## AGE_AT_DIAGNOSIS    0.036276  1.036942  0.002729 13.292   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                    exp(coef) exp(-coef) lower .95 upper .95
## ir_is                 0.8687     1.1512    0.6966     1.083
## Pam50_SUBTYPEBasal    1.1995     0.8337    0.9274     1.552
## Pam50_SUBTYPEHer2     1.3219     0.7565    1.0189     1.715
## Pam50_SUBTYPELumA     0.8051     1.2421    0.6431     1.008
## Pam50_SUBTYPELumB     1.1198     0.8930    0.8896     1.410
## Pam50_SUBTYPENC       1.1125     0.8988    0.4517     2.740
## AGE_AT_DIAGNOSIS      1.0369     0.9644    1.0314     1.043
## 
## Concordance= 0.611  (se = 0.009 )
## Likelihood ratio test= 234.7  on 7 df,   p=<2e-16
## Wald test            = 224.2  on 7 df,   p=<2e-16
## Score (logrank) test = 228.9  on 7 df,   p=<2e-16